#ifndef Algebra_sgPardiso_h
#define Algebra_sgPardiso_h

#include "AlgebraModule.h"
#include "LinearSolver.h"

SG_NAMESPACE_OPEN

namespace Algebra
{
    /**
     * @brief Wrap oneMKL PARDISO to solve sparse linear equations.
     *
     * The Intel oneAPI Math Kernel Library (oneMKL) PARDISO package is a high-performance,
     * robust, memory efficient, and easy to use software package for solving large sparse
     * linear systems of equations on shared memory multiprocessors.
     *
     * @see https://www.intel.com/content/www/us/en/docs/onemkl/developer-reference-c
     * /2024-2/onemkl-pardiso-parallel-direct-sparse-solver-iface.html
     *
     * @author nene
     * @date Augst, 2024.
     */
    // Extend for the complex data type
    template <typename TDataType>
    struct TPardisoPrivate;

    template <typename TDataType>
    class ALGEBRA_EXPORT TPardiso : public TLinearSolver<TDataType>
    {
      public:
        ///@name Type Definitions
        ///@{
        typedef /*typename*/ TDataType                 DataType;
        typedef /*typename*/ Vector_t<TDataType>       VectorType;
        typedef /*typename*/ Matrix_t<TDataType>       MatrixType;
        typedef /*typename*/ SparseMatrix_t<TDataType> SparseMatrixType;
        ///@}

        TPardiso ();
        ~TPardiso () override;

        const char* getClassName () override;

        /// perform a numeric factorization of the given sparse matrix.
        /// @note KEEP THE MATRIX ALIVE DURING THE WHOLE SOLVING PROCESS.
        void compute (const SparseMatrixType& matrix) override;

        /// evaluate the solution x of Ax=b using the current decomposition of A.
        void solve (const VectorType& b, VectorType& x) const override;

        /// Computes the inverse of square matrix `b` (must be invertible) and stores it in `x` using direct solvers
        void solve (const SparseMatrixType& b, SparseMatrixType& x) const override;

        /// report whether previous computation was successful.
        Info_t info () const override;

      private:
        /// reset the solver state
        void init ();
        /// de-allocate the resouces for the factorization.
        void clearFactors ();

        // void checkMatrix ();

      private:
        bool m_factorizationIsOk = false;
        /// take Adapter Pattern to hide the underlying details.
        struct TPardisoPrivate<DataType>* d = nullptr;
    };

    // // Enforce the explicit instantiations to hide the underlying implementation.
    // template class ALGEBRA_EXPORT TPardiso<Real_t>;
    // template class ALGEBRA_EXPORT TPardiso<Complex_t>;

    // // After that, export the routines which has been explicitly specialized.
    // // Otherwise, the compiler will complain that these routines could not be found.
    // template <>
    // ALGEBRA_EXPORT TPardiso<Real_t>::TPardiso ();
    // template <>
    // ALGEBRA_EXPORT TPardiso<Real_t>::~TPardiso ();
    // template <>
    // const char* ALGEBRA_EXPORT TPardiso<Real_t>::getClassName ();

    // template <>
    // ALGEBRA_EXPORT TPardiso<Complex_t>::TPardiso ();
    // template <>
    // ALGEBRA_EXPORT TPardiso<Complex_t>::~TPardiso ();
    // template <>
    // const char* ALGEBRA_EXPORT TPardiso<Complex_t>::getClassName ();

    // // alias for compatibilities with the external modules
    // typedef TPardiso<Real_t> sgPardiso;

}  // namespace Algebra

SG_NAMESPACE_CLOSE

#endif  // Algebra_sgPardiso_h